home *** CD-ROM | disk | FTP | other *** search
/ MacFormat 1994 August / August CD.bin / Shareware / Education / numericalmethods Folder / chap_2 / a2_1.m next >
Encoding:
Text File  |  1994-06-05  |  2.2 KB  |  101 lines  |  [MATS/MATL]
  1. echo off;
  2. % NUMERICAL METHODS: MATLAB Programs, (c) John H. Mathews 1994
  3. % To accompany the text:
  4. % NUMERICAL METHODS for Mathematics, Science and Engineering, 2nd Ed, 1992
  5. % Prentice Hall, Englewood Cliffs, New Jersey, 07632, U.S.A.
  6. % This free software is complements of the author.
  7.  
  8. % Algorithm 2.1  (Fixed Point Iteration).
  9. % Section    2.1,    Iteration for Solving  x = g(x), Page 51
  10. echo on; clc; format long; hold off; clear
  11. % This program implements fixed point iteration.
  12.  
  13. %    Define and store g(x)    in the M-file  g.m
  14. % function y = g(x)
  15. % y = 1 + x - x.^2 ./4;
  16.  
  17. delete g.m
  18. diary g.m; disp('function y = g(x)');...
  19.            disp('y = 1 + x - x.^2 ./4;');...
  20. diary off;
  21.  
  22. % Remark. g.m and fixpt.m are used for Algorithm 2.1
  23. g(0); % Test for file g.m
  24. pause    % Press any key to see the graph y = g(x).
  25.  
  26. clc; clg;
  27. a = 0;
  28. b = 5;
  29. c = 0;
  30. d = 5;
  31. h = (b-a)/150;
  32. X = a:h:b;
  33. Y = g(X);
  34. X1 = [-1000 1000];
  35. Y1 = [-1000 1000];
  36. axis([a b c d]);...
  37. plot(X1,Y1,'-r',X,Y,'-g');...
  38. hold on;...
  39. plot([a b],[0 0],'b',[0 0],[c d],'b');...
  40. xlabel('x');...
  41. ylabel('y');...
  42. title('The line y = x and the curve y = g(x).');...
  43. grid;...
  44. axis;...
  45. hold off;...
  46. shg; pause    % Press any key to continue.
  47.  
  48. clc;
  49. %    Place the starting value in  p0
  50.  
  51. %    Place the number of iterations in  max1
  52.  
  53. %    Place the tolerance in  delta
  54.  
  55. p0 = 4.0;
  56. max1  = 100;
  57. delta = 1e-9;
  58.  
  59. [pc,err,P] = fixpt('g',p0,delta,max1);
  60.  
  61. pause    % Press any key for the list of iterations.
  62.  
  63. clc; clg;
  64. max1 = length(P);
  65. for j = 1:max1-1,
  66.   k1 = 2*j-1;
  67.   k2 = 2*j;
  68.   Vx(k1) = P(j);
  69.   Vy(k1) = P(j);
  70.   Vx(k2) = P(j);
  71.   Vy(k2) = P(j+1);
  72. end
  73. Vy(1) = 0;
  74.  
  75. c = 0;
  76. d = 2;
  77. Z0 = zeros(1,length(P));
  78. axis([a b 0 2]);...
  79. plot(X1,Y1,'-g',X,Y,'-g',Vx,Vy,'-r',P,Z0,'or');...
  80. hold on;...
  81. plot([a b],[0 0],'b',[0 0],[c d],'b');...
  82. xlabel('x');...
  83. ylabel('y');...
  84. title('Graphical analysis for fixed point iteration.');...
  85. grid;...
  86. axis;...
  87. hold off;...
  88. shg; pause    % Press any key to continue.
  89.  
  90. max1 = length(P);
  91. J = 1:max1;
  92. points = [J;P];
  93. Mx1 = 'Computations for the fixed point iteration method.';
  94. Mx2 = '     k                  p(k)';
  95. Mx3 = 'The fixed point is g(p) = p = ';
  96. Mx4 = 'The error estimate for p is  ± ';
  97. clc,echo off,diary output,...
  98. disp(''),disp(Mx1),disp(''),disp(Mx2),disp(points'),...
  99. disp(''),disp(Mx3),disp(pc),...
  100. disp([Mx4,num2str(err)]),diary off,echo on
  101.